home *** CD-ROM | disk | FTP | other *** search
/ Best of www.BestZips.com (Collector's Edition) / Best of WWW.BESTZIPS.COM Collector's Edition (JCSM Shareware) (JCS Marketing).ISO / prgtools / tn2.zip / SIMPLEX.T < prev    next >
Text File  |  1996-11-15  |  5KB  |  320 lines
  1. %
  2. % "simplex.t" solves a linear programming problem 
  3. % using the simplex method
  4. %
  5. %   Sample program for the T Interpreter by:
  6. %
  7. %   Stephen R. Schmitt
  8. %   962 Depot Road
  9. %   Boxborough, MA 01719
  10. %
  11.  
  12. const N : int := 3      % number of variables
  13. const M : int := 4      % number of conditions
  14.  
  15. var Sx, Tx    : array[M+1,N+1] of real
  16. var basis     : array[M+1] of int
  17. var non_basis : array[N+1] of int
  18.  
  19. program
  20.  
  21.     var row, col : int
  22.     label program_exit :
  23.  
  24.     initialize
  25.  
  26.     loop
  27.  
  28.         show_array
  29.  
  30.         col := pivot_col
  31.         exit when col = 0   % all done
  32.  
  33.         row := pivot_row( col )
  34.         if row = 0 then
  35.  
  36.             put "infinite solution"
  37.             goto program_exit
  38.  
  39.         end if
  40.  
  41.         exchange( row, col )
  42.         update_array( row, col )
  43.  
  44.     end loop
  45.  
  46.     show_result
  47.  
  48.     program_exit:
  49.  
  50. end program
  51.  
  52. %
  53. % "initialize" sets up the simplex array for the constrained optimization
  54. % problem:
  55. %
  56. %   maximize:    y  =  2*x1 - 3*x2 + 3*x3
  57. %   subject to:  2 >= -2*x1 + 3*x2
  58. %                5 >=  1*x1 - 2*x2 - 4*x3
  59. %                6 >=  2*x1 + 1*x2 + 1*x3
  60. %               10 >=  1*x1 + 1*x2 + 5*x3
  61. %
  62. %   and all variables are non-negative
  63. %
  64. procedure initialize
  65.  
  66.     % variable indices
  67.     non_basis[0] := 0
  68.     non_basis[1] := 1
  69.     non_basis[2] := 2
  70.     non_basis[3] := 3
  71.  
  72.     % slack variable indices
  73.     basis[0] := 0
  74.     basis[1] := N + 1
  75.     basis[2] := N + 2
  76.     basis[3] := N + 3
  77.     basis[4] := N + 4
  78.  
  79.     % maximize y = 2*x1 - 3*x2 + 3*x3
  80.     Sx[0,0] :=  0
  81.     Sx[0,1] := -2
  82.     Sx[0,2] :=  3
  83.     Sx[0,3] := -3
  84.  
  85.     % 2 >= -2*x1 + 3*x2
  86.     Sx[1,0] :=  2
  87.     Sx[1,1] := -2
  88.     Sx[1,2] :=  3
  89.     Sx[1,3] :=  0
  90.  
  91.     % 5 >=  1*x1 - 2*x2 - 4*x3
  92.     Sx[2,0] :=  5
  93.     Sx[2,1] :=  1
  94.     Sx[2,2] := -2
  95.     Sx[2,3] := -4
  96.  
  97.     % 6 >=  2*x1 + 1*x2 + 1*x3
  98.     Sx[3,0] :=  6
  99.     Sx[3,1] :=  2
  100.     Sx[3,2] :=  1
  101.     Sx[3,3] :=  1
  102.  
  103.     % 10 >=  1*x1 + 1*x2 + 5*x3
  104.     Sx[4,0] := 10
  105.     Sx[4,1] :=  1
  106.     Sx[4,2] :=  1
  107.     Sx[4,3] :=  5
  108.  
  109. end procedure
  110.  
  111. %
  112. % "pivot_col" finds the column location for pivoting
  113. %
  114. function pivot_col : int
  115.  
  116.     var j, col : int
  117.     var t : real
  118.  
  119.     col := 0
  120.     t := 0
  121.  
  122.     for j := 1...N do
  123.  
  124.         continue when Sx[0,j] >= 0.0
  125.  
  126.         if t > Sx[0,j] then
  127.  
  128.             t := Sx[0,j]
  129.             col := j
  130.  
  131.         end if
  132.  
  133.     end for
  134.  
  135.     return col
  136.  
  137. end function
  138.  
  139. %
  140. % "pivot_row" finds the row location for pivoting
  141. %
  142. function pivot_row( col : int ) : int
  143.  
  144.     var i, row : int
  145.     var r, t : real
  146.  
  147.     row := 0
  148.     t := 1.0e+99
  149.  
  150.     for i := 1...M do
  151.  
  152.         continue when Sx[i,col] = 0.0
  153.         
  154.         r := Sx[i,0] / Sx[i,col]
  155.  
  156.         if r > 0 and r < t then
  157.  
  158.             t := r
  159.             row := i
  160.  
  161.         end if
  162.  
  163.     end for
  164.  
  165.     return row
  166.  
  167. end function
  168.  
  169. %
  170. % "exchange" switches the switches basis and non-basis variable indices
  171. % corresponding to the pivot row and column
  172. %
  173. procedure exchange( row, col : int )
  174.  
  175.     var t : int
  176.  
  177.     t := basis[row]
  178.     basis[row] := non_basis[col]
  179.     non_basis[col] := t
  180.  
  181. end procedure
  182.  
  183. %
  184. % "update_array" calculates the next set of values of the simplex array
  185. % for a given pivot row and column
  186. %
  187. procedure update_array( row, col : int )
  188.  
  189.     var pv, rv, cv : real
  190.     var i, j : int
  191.  
  192.     % update pivot value
  193.     pv := Sx[row,col]
  194.     Tx[row,col] := 1 / pv
  195.  
  196.     % update pivot row element values
  197.     for j := 0...N do
  198.  
  199.         continue when j = col
  200.         Tx[row,j] := +Sx[row,j] / pv
  201.  
  202.     end for
  203.  
  204.     % update pivot column element values
  205.     for i := 0...M do
  206.  
  207.         continue when i = row
  208.         Tx[i,col] := -Sx[i,col] / pv
  209.  
  210.     end for
  211.  
  212.     % update the rest of the simplex array
  213.     for i := 0...M do
  214.  
  215.         continue when i = row
  216.  
  217.         for j := 0...N do
  218.  
  219.             continue when j = col
  220.             Tx[i,j] := Sx[i,j] - Sx[row,j] * Sx[i,col] / pv
  221.  
  222.         end for
  223.  
  224.     end for
  225.  
  226.     % copy temporary array into simplex array
  227.     for i := 0...M do
  228.  
  229.         for j := 0...N do
  230.  
  231.             Sx[i,j] := Tx[i,j]
  232.  
  233.         end for
  234.  
  235.     end for
  236.  
  237. end procedure
  238.  
  239. %
  240. % "show_array" displays the simplex array
  241. %
  242. procedure show_array
  243.  
  244.     var i, j : int
  245.  
  246.     put "      "...
  247.  
  248.     for j := 0...N do
  249.  
  250.         put non_basis[j]:12...
  251.  
  252.     end for
  253.  
  254.     put ""
  255.  
  256.     put "      "...
  257.  
  258.     for j := 0...N do
  259.  
  260.         put "------------"...
  261.  
  262.     end for
  263.  
  264.     put ""
  265.     
  266.     for i := 0...M do
  267.  
  268.         put basis[i]:4, " |"...
  269.  
  270.         for j := 0...N do
  271.  
  272.             put Sx[i,j]:12...
  273.  
  274.         end for
  275.  
  276.         put ""
  277.  
  278.     end for
  279.  
  280.     put ""
  281.  
  282. end procedure
  283.  
  284. %
  285. % "show_result" shows the optimum value of the objective function
  286. % and the corresponding values of the variables
  287. %
  288. procedure show_result
  289.  
  290.     var i, j : int
  291.     var found : boolean 
  292.  
  293.     put "optimum y  = ", Sx[0,0], " at:"
  294.  
  295.     for j := 1...N do
  296.  
  297.         put "        x", j, " = "...
  298.  
  299.         found := false
  300.         
  301.         for i := 1...M do
  302.  
  303.             if basis[i] = j then
  304.  
  305.                 put Sx[i,0]
  306.                 found := true
  307.  
  308.             end if
  309.  
  310.         end for
  311.  
  312.         if not found then
  313.  
  314.             put 0
  315.  
  316.         end if
  317.  
  318.     end for
  319.     
  320. end procedure